home *** CD-ROM | disk | FTP | other *** search
/ Resource for Source: C/C++ / Resource for Source - C-C++.iso / codelib8 / v_10_08 / 1008050a < prev    next >
Encoding:
Text File  |  1995-11-01  |  3.0 KB  |  96 lines
  1. /*                    Listing 3                     */
  2. /*****************************************************
  3.                 Name: SRF_INCP.C
  4.          Description: Functions for calculating ray /
  5.                       surface intercepts
  6. Global Function List: srf_intrcpt
  7.          Portability: Standard C
  8. *****************************************************/
  9.  
  10. #include <stdlib.h>
  11. #include <float.h>
  12. #include <math.h>
  13. #include <srf_incp.h>
  14.  
  15. /*****************************************************
  16.        Name: srf_intrcpt
  17.  Parameters: Order - order of surface 1 (plane) or 2
  18.              PrvLoc - array of previous ray x,y,z loc.
  19.              CrtDir - array of current ray direction
  20.              Coef - matrix of surface coefficients
  21.              Disp - array of distances between
  22.              previous location and surface intercept
  23.      Return: Number of valid surface intercepts 0,1,2
  24. Description: Finds the surface ray intercepts are
  25.              found for either a first or second order
  26.              surface.  The distances to the surface
  27.              intercepts are stored in Disp and the
  28.              number of surface intercepts is returned.
  29. *****************************************************/
  30. int srf_intrcpt( int Order, double *PrvLoc,
  31.       double *CrtDir, double **Coef, double *Disp )
  32.    {
  33.    double a, b, c, SqrtArg = 0.0, Denom, Coefij;
  34.    size_t i, j;
  35.    Denom = Disp[0] = Disp[1] = DBL_MAX;
  36.    
  37.    /* Get the 0, 0 order terms. */
  38.    a = b = 0.0;
  39.    c = Coef[0][0];
  40.  
  41.    /* Get the 0, j order terms. */
  42.    for ( j = 0; j < 3; j++ )
  43.       {
  44.       Coefij = Coef[0][j + 1];
  45.       b += Coefij * CrtDir[j];
  46.       c += Coefij * PrvLoc[j];
  47.       }
  48.  
  49.    if ( Order == 2 )
  50.       {
  51.       /* Get the i, j order tems.  The Coef is out
  52.       ** of sync by 1 with the vector indcies. */
  53.       for ( i = 0; i < 3; i++ )
  54.          {
  55.          for ( j = i; j < 3; j++ )
  56.             {
  57.             Coefij = Coef[i + 1][j + 1];
  58.             a += Coefij * CrtDir[i] * CrtDir[j];
  59.             b += Coefij * ( PrvLoc[i] * CrtDir[j] +
  60.                   PrvLoc[j] * CrtDir[i] );
  61.             c += Coefij * PrvLoc[i] * PrvLoc[j];
  62.             }
  63.          }
  64.       SqrtArg = b * b - 4.0 * a * c;
  65.       Denom = 2.0 * a;
  66.       /* If the roots are imaginary return */
  67.       if ( SqrtArg < 0.0 ) return ( 0 );
  68.       }   /* if Order == 2 */
  69.  
  70.    /* If the denomiator a is very small treat the
  71.    ** surface as plane */
  72.    if (( fabs( a ) < DBL_EPSILON ) || ( Order == 1 ))
  73.       {
  74.       if ( fabs( b ) > DBL_EPSILON )
  75.          {
  76.          Disp[0] = -c / b;
  77.          return ( 1 );
  78.          }
  79.       return ( 0 );
  80.       }   /* if fabs( a ) */
  81.  
  82.    /* Special case one root. */
  83.    if ( SqrtArg < DBL_EPSILON )
  84.       {
  85.       Disp[0] = -b / Denom;
  86.       return ( 1 );
  87.       }   /* if SqrtArg */
  88.  
  89.    /* Normal case two real roots. */
  90.    SqrtArg = sqrt( SqrtArg );
  91.    Disp[0] = ( -b + SqrtArg ) / Denom;
  92.    Disp[1] = ( -b - SqrtArg ) / Denom;
  93.    return ( 2 );
  94.    }   /* funciton srf_intrcpt */
  95. /* End of File */
  96.